Epidemic spreading on interconnected networks 
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Many real networks are not isolated from each other but form networks of networks, often inter- 
related in non trivial ways. Here, we analyze an epidemic spreading process taking place on top of 
two interconnected complex networks. We develop a heterogeneous mean field approach that allows 
us to calculate the conditions for the emergence of an endemic state. Interestingly, a global endemic 
state may arise in the coupled system even though the epidemics is not able to propagate on each 
network separately, and even when the number of coupling connections is small. Our analytic results 
are successfully confronted against large-scale numerical simulations. 
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Epidemic spreading is one of the most successful ap- 
plication areas of the new science of networks [TJ [2] . In- 
deed, the general acceptance within the scientific com- 
munity that many diseases, like sexually transmitted dis- 
eases or the H1N1 virus, spread over networked systems 
represents a major step toward their understanding and 
control [3H5]. From a physics perspective, epidemic pro- 
cesses have been widely studied as a paradigm of non- 
equilibrium phase transitions with absorbing states [6]. 
When applied to complex networks, these processes have 
become a source of new and striking phenomena that 
do not have a counterpart in regular lattices. Germane 
examples are the absence of epidemic and percolation 
thresholds in scale-free networks with a power law de- 
gree distribution P(k) ~ fc -7 with 7 <G (2,3], and an 
anomalous critical behavior when 7 £ (3,4) [7HH]. 

We currently have a solid understanding of epidemic 
processes when they take place on single isolated net- 
works. In contrast, our comprehension is very limited 
when epidemics happen on coupled interconnected net- 
works. For example, sexually transmitted diseases can 
propagate both in heterosexual and homosexual networks 
of sexual contacts [3] . These two networks are not com- 
pletely isolated due to the existence of bisexual individ- 
uals, which act as an effective coupling between the two 
networks and potentially affect their epidemic proper- 
ties (TO]. To the best of our knowledge, a theory describ- 
ing these type of systems has not yet been fully devel- 
oped. 

In this paper, we fill this gap and present a rigor- 
ous heterogeneous mean field study of the susceptible- 
infected-susceptible (SIS) model taking place on two in- 
terconnected complex networks. Our analysis reveals a 
highly non-trivial behavior of the epidemic process de- 
pending on the strength and nature of the coupling be- 
tween the networks. We calculate the global epidemic 
threshold of the process, which turns out to be smaller 
than the epidemic thresholds of the two networks sep- 
arately under certain conditions. This implies that an 
endemic state may arise even if the epidemics is not en- 



demic in any of the two networks separately, as we prove 
analytically and with large-scale computer simulations. 

To begin our analysis, we have to specify the topologi- 
cal properties of our networks. Let A and B be two inter- 
connected random networks with given local properties 
and given two-points correlations. Each node in network 
A (respectively B) is characterized by a vector degree 
k a = {k aa ,k a b) representing the number of its internal 
connections k aa to other nodes in A, and the number of 
external connections k a b with nodes of network B. The 
analogous to the degree distribution in single networks is 
then the probability that a randomly chosen node of net- 
work A has a vector deg ree k a , P A (k a *j — PA(k aa: k a ]^), 
where, in general, k aa and k ab may be correlated. A sim- 
ilar definition applies to Pe(kb) for network B. These 
two distributions are arbitrary except for the consistency 
condition N A /N B = (k ba ) / (k ab ) , where N A /N B is the 
ratio between the sizes of the two networks and the right 
hand side is the ratio between their external average de- 
grees. This condition states that the total number of 
edges leaving network A towards B must be the same 
going from network B to network A. Two-points correla- 
tions, on the other hand, are encoded by the transition 
probabilities P AA (k' a \k a ), P AB (k' b \k a ), PB A (k^|k b ),and 
P BB (k' b \ki,). For instance, P^^lka) is the probability 
that, being in a node of network A with vector degree 
k Q a randomly chosen neighbor in network B has vector 
degree k' b . Similar definitions applies to the rest of the 
transition probabilities. 

As the epidemic spreading model, we consider the 
SIS model that, together with the susceptible-infected- 
recovered model (SIR), is one of the best studied models 
in epidemiology |llj . The model has a non-equilibrium 
phase transition between an endemic state with sustained 
epidemic activity and a healthy phase where the epi- 
demics dies out. Individuals can be in two different 
states, either susceptible (S) or infected (I). Infected in- 
dividuals decay spontaneously to the susceptible state at 
rate 5 (that, without loss of generality, we set to 8 = 1) 
whereas susceptible ones get infected at a rate propor- 
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tional to the number of infected neighbors they have at a 
given time. For two interconnected networks, we have to 
specify these processes separately. Let A aa (X bb ) be the 
infectious rate between nodes in network A (B) and A ab 
(A 6a ) the infectious rate from a node in A (B) to a node 
in B (A). 

The quantity of interest in the SIS dynamics is the 
prevalence, pit), defined as the fraction of infected nodes 
at a given time. To describe its evolution, we use the 
heterogeneous mean field approximation. In this approx- 
imation, nodes are classified within classes of equivalence 
such that all nodes within a given class are considered as 
statistically equivalent. In our case, classes of equivalence 
are defined by the network itself (A or B) and by specific 
values of the vector degree. Following these assumptions, 
we define the partial prevalences p£ (t) and as the 
fraction of infected nodes with a given vector degree. The 
total prevalences for each network is then computed as 
P A (t) = E ka PA(K)pt(t) and p B {t) = £ k6 P B (k 6 )p£. 
Following [51 H2] , the time evolution of the dynamics can 
be written as 



and 
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The first line in Eq. ([!]) describes the standard SIS model 
for a single network [12], whereas the second line ap- 
pears due to the coupling of network A with network B. 
An analogous equation can be written for network B by 
swapping the indices A — > B and a — >■ b in Eq. . 

As in a single network, this process undergoes a phase 
transition between a healthy phase with p£ = p^ = 
and an endemic phase with pj^ ^ and p^ b ^ 0. How- 
ever, a mixed phase with endemic activity in one network 
whereas the other is in a heathy state is not possible in 
the system formed by the two coupled networks, where 
the epidemics propagates to the whole system if it is able 
to propagate in one of the networks. This is due to the 
fact that the state p£ = and p^ b ^ is not a fixed 
point of the dynamics in Eq. . The critical point sep- 
arating the healthy and endemic phases can be obtained 
by studying the stability of the absorbing solution. This 
can be done by linearizing the system of Eqs. ([I]) around 
p A = Pk b = and studying the spectrum of the corre- 
sponding matrix. Close to the absorbing state, Eqs.Q 
can be written as 



dp 
di 



(2) 



\ aa k aa p AA (K\k a ) x ba k ab p AB (k' b \k a 



\ ab h a p BA (k' a \k b ) \ bb k bb p BB (k' b \k b ) 



(3) 



The absorbing state is stable whenever the maximum 
eigenvalue of matrix C satisfies A m < 1. Otherwise, the 
absorbing state is unstable and an endemic state takes 
over the system. The critical epidemic point is then de- 
fined by A TO = 1. 

When the networks do not have two-points correla- 
tions, the transition probabilities can be written as 



Paa(KK) 



fcaa-P(ka) p ,,,,, ^ k' ba P(k' b ) { A\ 
JjT^ ' P-AB(k fc |k Q ) — jj—^ (4) 



In this case, the eigenvalues of matrix C are the solutions 
of the following equation 
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and where 



where we define the vector prevalence as p = (p^ , p£ ) 



_ \aa\ab (k aa k ab ) _ . bb . ba (k bb k ba ) 

[kaa) [kab) {fob) [ha) 

Constants appearing on the left hand side of Eq. ^ have 
a clear interpretation. Indeed, A A and A B are the max- 
imum eigenvalues of networks A and B as if they were 
isolated. Therefore if, for instance, A A > 1 then network 
A is able to sustain an endemic state by itself when iso- 
lated from network B. Similarly, A AB is the maximum 
eigenvalue of the network AB as a pure bipartite system, 
that is, when all internal connections inside networks A 
and B are absent. Again, when A AB > 1, the pure bipar- 
tite network AB is able to sustain an endemic state, both 
in the coupled system and even if there were no connec- 
tions whatsoever within each individual network. Con- 
stants a ab , a ba , Pab and p ba appearing on the right hand 
side of Eq. |5]) contain information about the strength 
and nature of the coupling between the nets. 

From Eq. ([5]), it is easy to see that the maximum 
eigenvalue of matrix C, A m , is always larger than 
max (A^, As, A^s). It is therefore possible to find en- 
demic states with A m > 1 but where A^ < 1, A^ < 1, 
and A AB < 1, that is, situations where neither networks 
A and B isolated nor the pure bipartite network AB are 
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able to sustain the endemic state and yet the epidemics 
pervades in the coupled system. This effect is more or 
less important depending on the strength of the coupling, 
i. e., the number of links between the two networks - 
quantified by (k ab ) and (k ba )- and the specific correla- 
tions between internal and external degrees, measured 
by the factors (k aa k ab ) and (k ba k bb ). 

In the rest of the paper, we focus on the interesting case 
A a < 1, Ab < 1, and Aab < 1, and ask under which 
condition the endemic state exists, that is , A m > 1. 
From Eq. (f5j), we see that this happens when the right 
hand side of Eq. ^ evaluated at x = 1 is larger than the 
left hand side evaluated at the same point, which after 
some algebra yields 

[a ab + [i ab (l-A A )}[a ba + n ba (l-A B )} > (1- A A )(1- A B ). 

(9) 

This equation is one of the main results of our paper. It 
allows us to evaluate the conditions for the emergence of 
the endemic state in many different situations. 

In real networks of sexual contacts, the most promis- 
cuous individuals in one network are also the ones with 
the largest number of sexual partners in the other net- 
work [10 . This fact suggests that there exists a posi- 
tive correlation between the internal and external degree 
for a given node. We model these correlations by tak- 
ing the vector degree distribution to be PA(k a a,k ab ) — 
PA{kaa)g(k a b\kaa) , _where g(k ab \k aa ) is a Poisson distri- 
bution with mean k ab (k aa ) — (k a b)k^ a /{k^ a ), and analo- 
gously for PB{k bb , k ba ). This choice allows us to interpo- 
late between a random distribution of links between the 
two networks, (3 = 0, and positive correlations, where 
high degree nodes in both networks concentrate the ma- 
jority of the coupling links, (3 > 0. Inserting this assump- 
tion in Eq. ([9| , we obtain the critical lines that define the 
phase diagram in the hyperplane A^, Ab, ot ab , and a ba . 

Figure [T] shows examples of phase diagrams for two 
identical networks with a symmetric coupling, that is, 
Aa = A B = A, a ab = a ba = a, and fi a b = Mba = M- In 
this case, Eq. ^ simplifies to a > (1 — A)(l — /i). We 
consider the case of networks with an exponential degree 
distribution and the values of (3 = 0, 1, 2. The most inter- 
esting case corresponds to the range of parameters where 
both internal and external links are needed for the ex- 
istence of the endemic state -depicted in orange in Fig. 
[T] and particularly when the internal connections out- 
number the external ones (the portion below the dashed 
line). This area shrinks when increasing but, simul- 
taneously, it also appears at lower values of the ratio 
of external versus internal connections for the same A. 
Therefore, keeping the two networks unchanged and for 
a fixed number of links between them, the epidemic can 
either be in the healthy or endemic phase depending on 
how these links are distributed among the nodes of the 
networks. Notice that, if correlations are strong enough, 
the bipartite network alone is able to sustain by itself the 
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FIG. 1: Phase diagram showing the healthy phase and the 
endemic phase for the case of a random coupling (top), linear 
correlations (middle) and super-linear correlations between 
external and internal degrees (bottom). In all cases, networks 
A and B are identical, with an exponential degree distribu- 
tion with (k aa ) — (kbb) = 10 and minimum degree 2. The 
area depicted in gray corresponds to the healthy phase and 
the transition to the endemic state is represented by the black 
curve. The orange area corresponds to the endemic state that 
needs the contribution of both internal and external links and 
the white area to the case when the pure bipartite network 
alone is able to sustain the endemic state. The horizontal 
dashed line marks the point where external connections out- 
number internal ones. 

endemic state in the coupled system. 

We checked our predictions with large scale numerical 
simulations. The SIS dynamics is simulated with a con- 
tinuous time dynamics as follows. During the course of 
the simulation, we keep track of the number of infected 
nodes Nj(t) and the number of active links EA{t), where 
an active link is defined as a link connecting a suscepti- 
ble and an infected node. At each step, with probability 
p r = Ni(t)(Ni(t) + A-?a(*)) _1 ) a randomly chosen in- 
fected node is turned susceptible whereas, with probabil- 
ity 1 — p r , an active link is chosen at random and the sus- 
ceptible node attached to it is turned infected. After this 
procedure, time is updated as t — > t+(N](t) + XEA(t))^ 1 . 
We run this algorithm on the networks used in Fig. [I] of 
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FIG. 2: Time evolution of the prevalence p(t) below and above 
the critical point for A = 0.8 for two identical networks inter- 
connected with linear correlations, i. e., f3 — 1. Each network 
has an exponential degree distribution with (k) = 10 and 
minimum degree 2. The size of each network is N = 10 6 
and results are averaged over 100 different realizations of the 
process. Right at the critical point, the prevalence decays as 




FIG. 3: Critical lines for f3 = 0, 1, 2 for the same networks as 
in Fig. [2] Dotted lines correspond to numerical simulations, 
black lines are the theoretical predictions according to Eq. (J9|, 
whereas blue dashed lines are obtained by rescaling the A axis 
as A^A(l-(fc)/(fc 2 ))- 1 . 

size N = 10 6 . Having fixed the internal network prop- 
erties, we fix the value of A by adjusting the infectious 
rate A using Eq. Q. For each value of A, we change the 
average external degree and study the temporal behavior 
of the prevalence. The critical point is determined as the 
point where p(t) decays as a power law [6], as shown in 
Fig.|2] 

Using this method, we compute the critical line in 
the plane ((k a b) / (k aa ) , A) for the three cases analyzed 
in Fig. [I] We limit our search to the domain ((k a b) < 
(k aa ),A < 1). Simulation results are shown in Fig. [3] 
along with the analytic predictions given by Eq. We 
observe a systematic shift between theory and simula- 



tions. The reason is that our mean field approach does 
not consider dynamical correlations. These correlations 
imply that, with high probability, an infected node has 
its infecting neighbor still infected during some time right 
after the infectious event. This reduces the number of 
potential new infections that this node can produce. We 
can correct this effect in the same way as it is done in the 
SIR model, just by reducing the number of contacts by 1. 
In our equations, this is achieved by replacing the maxi- 
mum eigenvalue of the networks by A = X(k(k — 1)) /(k). 
This is, of course, a limiting case because there are cases 
where the infecting neighbor recovers and can be rein- 
fected. Thus, we expect to find simulations results be- 
tween these two extremes. We show this correction in 
Fig. [9] as the gray area. Indeed, all simulation points fall 
within this gray area, confirming then our intuition. 

The study of interconnected and/or interdependent 
networks reveals new and unexpected phenomena |13j . 
Here, we have shown that two networks well below their 
respective epidemic thresholds may sustain an endemic 
state when coupling connections are added, even in small 
number. This may have important implications for the 
design of efficient control strategies. However, the ef- 
fects of the coupling are highly non-trivial and may vary 
depending on the strength and the correlations of the 
interconnecting links. We foresee similar effects appear- 
ing in many different dynamics showing equilibrium and 
non-equilibrium phase transitions. 
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Note added.-During the final writing of this paper, we 
became aware of a recent preprint where the SIR model 
on interconnected networks is studied 14 1. In that work, 
the authors find a mixed phase where one network prop- 
agates the epidemic while the other does not. Neither 
our analytic nor simulation results indicate the existence 
of such mixing phase in our case. 
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